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1 Introduction 

There has been considerable interest in the use of 
Bayesian formulations and approaches in the 
solution of the biomagnetic inverse problem 
[1,2,3,4], The Bayesian Power Imaging (BPI) 
method is a new Bayesian approach that is 
introduced in [5]. In this companion paper, the BPI 
method is extended to analyze not one but two sets 
of experimental data in order to highlight the 
differences in their underlying source structures. 

2 Methods 

The BPI method as described in [5] is a fusion of the 
Bayesian approach detailed in [4] and the direct 
power algorithm derived in [6], The method 
provides a means of calculating the probability that 
features of an activation map are due to meaningful 
sources as opposed to noise. 

In order to highlight aspects of the BPI method that 
are central to the comparison calculations set out in 
this paper, we will summarise the relevant features 
of [5], The first step in the BPI method is to 
determine the a-posteriori distribution of the 
expected measurements for each channel, given the 
a-priori probability distribution of the current and 
the data. 

The means of this distribution are then used in the 
direct power method of [6] to calculate the 
activation A of the 0 th voxel; 

A e( m )=T J ( mT Ve) 

where I \r\ =(p 2 +d) PL k e and L k e is the 
measurement vector due to a current dipole at the 
position labelled 0 and direction k. 

P and D are the Gram-Schmidt matrix and the 
covariance matrix respectively. 

In order to determine the probability distribution on 
the activation we use the a-posteriori distribution on 
the measurement cha nn els. This results in quadratic 
forms in Gaussian random variables which are 
evaluated using Monte-Carlo integration. 


The calculation proceeds on a voxel by voxel basis 
and generates a significance or probability map by 
comparing the activation distribution due to the data 
with that due to random noise. The value at each 
voxel of the map is the probability that the activation 
at the voxel cannot be explained (to a chosen 
statistical level) by the noise. 

This probability is given by; 

p e ( A e Ue (m)-T - A e (z, ) ] 

Here A e (m) and A e [z.) are the activations of 
voxel 6 due to the measurement vector m and the 
i th noise vector z.. <t> denotes the Heaviside step 
function and N is the number of trials or noise 
vectors in the ensemble. A threshold T(A g (m)) is 
included to allow the exclusion of large regions of 
low significance from the resulting map. 

The probability P e is a measure of the significance of 
the activation on a voxel compared with that 
produced by noise. In effect, we are comparing the 
data against a prior model which, in this instance, is 
characterized by random noise. In so far as random 
noise simulates uncorrelated detector noise this prior 
model could be considered as equivalent to an 
idealised ‘empty room’ experiment. 

Other prior models can be used. Data from real 
'empty room' experiments would allow the 
examination of the significance of a response within 
the context of environmental noise and instrumental 
artifacts. Similarly, using data from a pre-stimulus 
period, would help to identify significance to 
features corresponding to the evoked response 
within the context of other subject derived sources, 
both cerebral and non-cerebral. 

This paper is focused on the use of the BPI to 
determine the significance of differences in the 
responses recorded in two otherwise identical 
experiments involving different stimuli or tasks. 
This is a very common paradigm in brain research 
where a variation of stimulus is used to identify a 
specific aspect of brain function. Using the BPI 
formulation, one experiment becomes the prior 
model against which the other is compared. The map 



of significant differences then indicates areas where 
the prior model is not a good explanation of the 
experiment, i.e. where the responses are different. 
There is a complication. The Monte-Carlo 
calculation requires the comparison of the data 
against many datasets representative of the prior 
model. Here our prior model consists of a single 
dataset. In order to proceed we subtract the 
activation of the prior model from that of the 
experiment to give a difference map. We then 
evaluate the significance of this difference map 
against the random noise as previously to determine 
the probability that the differences are not due solely 
to noise. If the difference map is denoted by 

A A e = A g ( m 1 ) -A e ( m 2 ) 

then 

Pg (A Ag ) = ^L^lA^ -T~Ag (*, )] 

Either of the two one-sided probabilities can be 
found by switching the role of prior model and 
experiment between the two datasets. 

3 Results 

This extension of the BPI method is illustrated by 
real and simulated data studies. For both studies the 
experimental instrument was the Neuromag 122™ 
full head detector system which consists of 61 pairs 
of first order gradiometers that measure the 
tangential gradients of the magnetic field. 

3.1 Simulated data study 

The simulated data was generated by two current 
dipoles embedded in a homogeneous sphere. The 
dipole magnitudes were the same. The source space 
was a part spherical shell surface of radius 8 cm 
which covered a 2 radian by 2 radian solid angle in 
posterior regions of the brain. The surface was 
divided into 33 by 33 voxels. The arrangement is 
shown in Figure 1. 



Figure 1 The position of the source space relative to 
the detector helmet. 

The first dataset m l was generated with the dipoles 
at positions (5cm, 12cm) and (10cm, 4cm) relative 
to the bottom left hand comer of the source space 
and oriented along the x and y axes respectively. For 
the second dataset, m 2 , the dipole at (5cm, 12cm) 
was moved 3 cm along the x axis to (8cm, 12cm) 
without any change in orientation. Uncorrelated 
Gaussian noise was added to both datasets so that 
the resulting signal to noise ratio was 10.2. 

As in [5], the regularisation parameter g was 
calculated as O.lxtrac e(P)/M where P is the 
Gram-Schmidt matrix and M is the number of 
detectors. The number of trials N in the Monte- 
Carlo calculation was 1000. The threshold T was set 
to two standard deviation above the average 
activation over the source space. 

Figure 2 shows plots, over the full source space, of 
the activations A g (m t ) and A g {m 2 ) and the two 
one-sided probability estimates of the differences. 

3.1 Real data study 

The real data chosen for the second illustration came 
from an experiment designed to elicit differences 
between the responses to different stimulus classes 
[7], Subjects were shown pictures of human faces, 
animals, motor bikes and abstract patterns in a 
random sequence. The task was to identify 
membership of a stimulus class by a cued key press 
or deny membership by pressing an alternative key. 



Figure 2 The two upper maps are the activation 
maps A g {m i ) and A g ( m 2 ) corresponding to the 
initial (left) and final (right) positions of the dipoles. 
The left lower map shows the map of significant 
differences derived from one sided probability 
calculations. For the left map the first data set was 
compared with a prior defined by the second data 
set. The right lower map was generated by reversing 
the roles of data and model. Each map is a two- 
dimensional projection of the source space. The 
dipole positions are indicated by asterisks. 

The goal was to identify any response that was face 
specific. We have reported such a face specific 
response at ~150 ms after stimulus onset, with the 
main centre of activation in inferior occipito-temoral 
cortex [7], For the present study we used averaged 
data from one subject who demonstrated such face 
specificity weakly. MRI data was available for this 
subject. 

The BPI method was applied using a source space 
consisting of the entire brain. The activation maps 
for the latency corresponding to the peak in the 
global signal power (160 ms) are complex with 
significant activity widespread over the occipital, 
temporal and parietal regions. The difference 
algorithm was also used with the responses to 
animals in the role of the prior model. 

Figure 3 is derived from the one-sided probability 
that the activation of the face response is 
significantly greater than the activation of the animal 
response. For simplicity the figure shows those 
regions where, to a high probability, the animal data 
prior does not explain the face response. 



Figure 3 The regions on the brain surface where the 
response to the face stimuli cannot be explained by 
the prior animal response model. 

4 Discussion 

The simulation study provides a simple illustration 
of the use of the BPI method to identify activation 
differences. The activation maps correspond closely 
to the known dipolar source structure with the usual 
smearing introduced by the choice of basis and 
regularisation. The thresholded maps of significant 
differences identify correctly the regions where 
there are significant differences between the two 
experiments. It should be noted that the source 
regions linked to the stationary dipole do not appear 
in the map of significant differences. 

The illustration based on the real data has provided a 
means of quantifying face specificity in source 
space. Figure 3 shows that, for this subject, face 
specificity is associated with right occipito-temporal 
regions. Both inferior and superior regions are 
indicated. This is broadly in line with previous 
source analyses which were reported alongside a 
statistical analysis of face specificity carried out in 
data space [7]. 

The identification of the anatomical regions in 
Figure 3 is not robust. The full analysis of the 
sources will require further refinement of the source 
space as well as careful optimisation of the 
regularisation conditions, the weighting distribution, 
and the noise covariance assumptions. However, the 
example does illustrate how the BPI method can be 
used to identify source regions where the activity is 
a significant function of the stimulus. It is worth 
noting that, once again, the strong occipital activity 
found in both the face and animal responses 
disappear in the map of the significant differences. 





There is considerable interest in the comparison of 
MEG images with the results obtained by other 
functional imaging methods, for example fMRI. The 
BPI method is used here in a directly analogous way 
to the subtraction techniques employed in fMRI. 
This enables easier comparison between the 
modalities. 
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